home *** CD-ROM | disk | FTP | other *** search
/ Pascal Super Library / Pascal Super Library (CW International)(1997).bin / MATH / NRPAS13 / CONVLV.PAS < prev    next >
Pascal/Delphi Source File  |  1991-04-29  |  2KB  |  46 lines
  1. PROCEDURE convlv(data: glnarray; n: integer; respns: glnarray; m: integer;
  2.        isign: integer; VAR ans: gln2array);
  3. (* Programs using routine CONVLV must define the types
  4. TYPE
  5.    glnarray = ARRAY [1..n] OF real;
  6.    gln2array = ARRAY [1..n2] OF real;
  7. where n is the dimension of the data and n2=2*n. NOTE: when used with CONVLV,
  8. the data dimension in FOUR1 and in TWOFFT must be the same as gln2array here.
  9. i.e. TYPE  gldarray = gln2array; gl2narray = gln2array *)
  10. VAR
  11.    no2,i,ii: integer;
  12.    dum,mag2: real;
  13.    fft: gln2array;
  14. BEGIN
  15.    FOR i := 1 TO ((m-1) DIV 2) DO BEGIN
  16.       respns[n+1-i] := respns[m+1-i]
  17.    END;
  18.    FOR i := ((m+3) DIV 2) TO (n-((m-1) DIV 2)) DO BEGIN
  19.       respns[i] := 0.0
  20.    END;
  21.    twofft(data,respns,fft,ans,n);
  22.    no2 := n DIV 2;
  23.    FOR i := 1 TO (no2+1) DO BEGIN
  24.       ii := 2*i;
  25.       IF (isign = 1) THEN BEGIN
  26.          dum := ans[ii-1];
  27.          ans[ii-1] := (fft[ii-1]*ans[ii-1]-fft[ii]*ans[ii])/no2;
  28.          ans[ii] := (fft[ii]*dum+fft[ii-1]*ans[ii])/no2
  29.       END ELSE IF (isign = -1) THEN BEGIN
  30.          IF ((sqr(ans[ii-1])+sqr(ans[ii])) = 0.0) THEN BEGIN
  31.             writeln('pause in routine CONVLV');
  32.             writeln('deconvolving at response zero'); readln
  33.          END;
  34.          dum := ans[ii-1];
  35.          mag2 := sqr(ans[ii-1])+sqr(ans[ii]);
  36.          ans[ii-1] := (fft[ii-1]*ans[ii-1]+fft[ii]*ans[ii])/mag2/no2;
  37.          ans[ii] := (fft[ii]*dum-fft[ii-1]*ans[ii])/mag2/no2
  38.       END ELSE BEGIN
  39.          writeln('pause in routine CONVLV');
  40.          writeln('no meaning for ISIGN'); readln
  41.       END
  42.    END;
  43.    ans[2] := ans[n+1];
  44.    realft(ans,no2,-1)
  45. END;
  46.